3D Dune Skeleton Model as a Coupled Dynamical System of 2D Cross-Sections 
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To analyze theoretically the stability of the shape and the migration process of transverse dunes 
and barchans, we propose a skeleton model of 3D dunes described with coupled dynamics of 2D 
cross-sections. First, 2D cross-sections of a 3D dune parallel to the wind direction are extracted 
as elements of a skeleton of the 3D dune, hence, the dynamics of each and interaction between 
them is considered. This model simply describes the essential dynamics of 3D dunes as a system of 
coupled ordinary differential equations. Using the model we study the stability of the shape of 3D 
transversal dunes and their deformation to barchans depending on the amount of available sand in 
the dune field, sand flow in parallel and perpendicular to wind direction. 
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Sand dunes, which are the largest granular objects on the Earth, move by wind and exhibit various morphodynamics. 
As typical shapes of dunes, barchan, transverse dune, linear dune, star dune, dome dune and parabolic dune are 
known [l], [^. The steadiness of wind direction and the amount of available sand in each dune field are considered as 
dominant factors determining these shapes. For example, unidirectional steady wind generates barchans or transverse 
dunes. The former are crescentic shaped dunes, and are formed in dune fields with small amounts of available sand, 
whereas transverse dunes which extend perpendicular to the wind direction, are formed in dune fields with the larger 
amount of available sand than the barchan-rich field. A characteristic aspect of recent dune studies is that quantitative 
analysis of dune morphodynamics has largely progressed. In particular, water tank experiments and computer models 
have uncovered the complex processes of dunes [3|-|8|. They were successful to reproduce formation and migration 
processes of barchans and other types of dunes under controlled setups. However, theoretical methodology to explain 
the complex morphodynamics of dunes beyond only reproducing them is yet to be developed. We, here, propose a 
skeleton model of 3D dunes described with coupled dynamics of 2D cross-sections j9|, which has a form of coupled 
ordinary differential equations. Using the model we study the morphodymamics of dunes, particularly the stability 
of the shape of transverse dunes and its deformation to barchans. 

Recently Katsuki and Nishimori has proposed a model for the collision dynamics of two 3D barchans focusing on 
the dynamics of their central 2D cross sections [8] [9]. That is called ABCDE(Aeolian/Aqueous Barchans Collision Dy- 
namical Equations). To introduce the present skeleton model for 3D transverse dunes, we employ similar assumptions 
to those used in ABCDE; 
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FIG. 1: (a) Outlook of skeleton model; triangular 2D-CSs (1 < i < A'^) constitute a skeletonized 3D transverse dune or a 3D 
barchan. (b) Shape similarity of 2D-CSs; independent of their size. Slope angles of d and ip are kept constant, then, horizontal 
position and the height of ith 2D-CS are uniquely determined if the coordinate of its crest {xi, hi) is given. 



First, laterally arraying 2D wind directional cross-sections (hereafter 2D-CSs) of a 3D transverse dune (or a barchan) 
are set as elements of the present skeleton model. Hence, the dynamics of each 2D-CS and the interaction between 
them is considered (Fig. [D^a)) 

As mentioned above, 3D barchans are observed in dune fields with small amount of available sand, thus, they are 
isolated on a hard ground both in wind direction and lateral direction. On the other hand, 3D transverse dunes 
observed under the same wind condition, extend in lateral direction therefore are not isolated. However, in wind 
direction, successive two laterally extending crests of 3D transverse dunes are separated by inter-dunes hard ground 
if the amount of available sand in the dune field is insufficient to cover all the ground, and we treat such cases 
here. Shortly, we treat 2D-CSs constituting 3D barchans or 3D transverse dunes which are isolated in wind direction 
from windward and leeward 2D-CSs and each of them is assumed to have triangular shape located on hard and flat 
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(a) intra-2D-CS (b) inter-2D-CS 




FIG. 2: (a) Intra-2D-CS sand flow; Over-crest sand flux q and sand trapping efficiency Tsihi) govern the intra-2D-CS sand 
flow, (b) Inter-2D-CS sand flow; The dominant parts of the flow in the upwind and downwind slopes are proportional to the 
height difference of neighboring slopes multiplied by the overlapped length of these slopes as indicated by grey areas. The 
additional contributions to the flow at the foot regions are indicated by black areas. 



ground. In addition, considering the fact that 2D-CSs of transverse dunes and barchans roughly have shape similarity 
independent of their size, we assume that the angles of their upwind and downwind slopes {9 and ip, respectively) are 
constant (Fig. [Hb)). Then, geometrical constants A, B, C are introduced like, 

tan 6* tan 09 „ tan 03 ^ tan 6* , , 

tan 9 + tan ip ' tan 9 + tan ip ' tan 9 + tan cp ' 

where, the values of A, B and C are set j^, | and i respectively reflecting typical 2D-CS profiles of real barchans and 
transverse dunes. Based on above assumptions, the horizontal (i.e., wind directional) position and the height of each 
2D-CS are uniquely determined if the coordinate {xi, hi) {1 < i < N) of its crest is given(Fig. [Hb)). 

Now, we consider the interaction between laterally neighboring 2D-CSs and then after the migration process of 
each. Dune migration occurs by the sand flow along the surface. In our description the sand flow is divided into two 
types: (a) intra-2D-CS flow (b) inter-2D-CS flow (Fig. IHa)). In the intra-2D-CS flow, the most relevant quantity for 
the present modelling is the over-crest sand flux q of each 2D-CS, which governs the erosion rate of its upwind slope. 
The over-crest sand flux q also determines the deposition rate of the downwind slope, and it is noted that finite ratio 
1 — Te (0 < 1 — Te < 1) of g is assumed to directly escape from the 2D-CS to the leeward inter-dune ground and 
the remaining ratio is deposited in the downwind slope. (Fig. Ufa)). This ratio Te is termed as the sand trapping 
efficiency, according to Momiji and Warren[10] . 

The over-crest flux q is assumed constant independent of the height of the 2D-CS, while Te is an increasing function 
of the height with limits lim Tsih) = 0, lim TE{h) ~ 1, with the specific form, 

TE{h) = h/{a + h),a^ 1.5. 

Next, the inter-2D-CS flux Du{i^j) / D^K^i^j^ between upwind/downwind slopes of neighboring 2D-CSs, i and j, is 
given as a function of heights and wind directional positions of their crests; 

According to previous 3D dunes formation modelsfllj. we approximate that these inter-2D-CS fluxes are primarily 
contributed by the lateral diffusion fluxes integrated over upwind and downwind slopes respectively, thus, are propor- 
tional to the height difference between neighboring slopes multiplied by the overlapped length of these slopes. These 
quantities are shown as grey parallelograms in Fig. Hfb). 

In addition, we consider the contributions of the lateral diffusion at the foots of the upwind and downwind slopes 
where the diffusion flux is not proportional to the height difference between neighboring slopes but to the local height 
of the corresponding slopes because of no neighboring 2D-CS exists at these regions. The diffused sand from these 
regions is assumed to be soon absorbed in the slopes of neighboring 2D-CS in short time. These quantities are shown 
as black triangles in Fig. [Hb). 

Note that we treat, in the below, only the cases of small height difference between neighboring slopes, therefore, 
the additional rule of diffusion at the foots-region of 2D-CSs docs not affect the following results. 

Specific forms of D^^^^j^ and D^(^j__^j^ are, 

where and Vd are phenomenological parameters to control the amount of inter-2D-CS sand flow at respective sides 
of slopes. 
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Now we consider the migration of each 2D-CS. As mentioned above, the size and the wind directional position of 
zth 2D-CS is uniquely determined if the coordinate {hi,Xi) of its crest is given. Therefore to describe the dynamics 
of {hi, Xi) (1 < i < N) corresponds to give the skeletonized dynamics of a 3D dune. 

Here, Aa;„i and lS.Xdi are, respectively, the wind directional displacement of the upwind and the downwind slopes 
of ith 2D-CS during Ai (Fig. [3ja)), while the change of the coordinate {hi,Xi) of the 2D-CS crest within the same 
interval are denoted as Ahi and Aa;^, respectively. Then, Ahi and Axi are expressed by Axdi and Axu 



Ah., = kAxdi - kAxui (3a) 
Ax, = BAxdi + GAxui (3b) 

where A, B, C are the geometrical constants introduced in ([1]). Consequently the eroded area per At along the upwind 
slope on ith 2D-CS is expressed as, 

2K + A/t» ^^^^ ^ ^^^^ ^ + + fD (4) 

where q is the over-crest sand flux, Z3„(i_yi_i) and Du(i^i+i) are inter-2D-CS flux from ith 2D-CS upwind slope to 
{i — l)th and {i + l)th 2D-CSs, respectively. The quantity /™ is incoming flux to ith 2D-CS from the windward 
intcr-dune ground. Similarly the deposited area per At along the downwind slope on ith 2D-CS is expressed as, 

2hi + Ahi ^ At{qTE{hi) + + Da(i^i+i)) (5) 

where Tsihi) is the sand trapping efficiency, Dd{i^i-i) and Dd(i^i+i) are the 2D-CS flux from ith 2D-CS downwind 

slope to (i— l)th and (i + l)th 2D-CSs, respectively. Using eqs. ((3a|)-([5|) and taking their limits At — > and Ahj — > 0, 
a system of coupled ordinary equations, 

^ = A L^TEihi) -1)+ J2 (^d(i-i) - Ou(i^J)) + /r ) (6a) 



' \ 3 = »±1 



+ cD„(,_j))-/rj (6b) 

(1 < j < N). 



is obtained, which describes the dynamics of a skeletonized 3D dune. 

Because dunes treated here are assumed to migrate on a flat and hard ground, its minimum height is zero that is 
not cared in (|6al) .(j6b p . Therefore we add a rule for the vanishment of 2D-CSs, that is, if hi decreases to hi = 0, then 
ith 2D-CS is taken as vanished. Note that the wind directional position x of already vanished 2D-CS is virtually kept 
allocated at the foot of the downwind foot of the neighboring 2D-CS in order to determine the sand flux from the 
neighboring 2D-CSs according to eqs. ([Sa)) and (|2b| (Fig. El^b)). 

Numerical simulation of eqs. (pa|) . (j6bl) is conducted with N = 1000 of 2D-CSs. The lateral boundary condition is 
set periodic, whereas the boundary condition in wind direction is set more carefully. Because the wind directional 
position X is variable of (|6bp thus boundary condition for the 2D-CSs is not required. However, as mentioned above, 
finite ratio 1 — Tsihi) of the over-crest flux q escapes from each 2D-CS. Then, to keep the realistic situation of dunes 
dynamics in desert fields where incoming sand flux to a dune is supplied by escaping sand flux from windward dunes. 




FIG. 3: (a) Combinatorial action of intra-2D-CS flow q and inter-2D-CS flow causes the erosion of upwind slope with width 
Axui and the deposition of downwind slope with the width Axdi- (b) Even after ith 2D-CS is vanished, its crest position Xi is 
virtually kept allocated to define the inter-2D-CS flow from neighboring 2D-CS. 
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the total amount of escaping sand, 

N 
1=1 

is set uniformly redistributed to the upwind slopes of each 2D-CS. It means that /™ — Ftotai/N in ((6a|) and (j6bp . 
Note that if finite number of 2D-CSs have already vanished they do not catch the redistributed sand, meaning that 
total amount of sand constituting 2D-CSs in the system is conserved only if all 2D-CSs are kept non-vanished, in that 
case, total area of 2D-CSs 

is kept constant too. 

In the present model of eqs. (|5a|) . (j6bp accompanied with eqs. (Pa|) . (j2bp . three-environmental parameters v^^Vd^q 
are introduced all of which increase if the wind force increases. Among them, Vd are related to the increase of sand 
flow in the lateral direction, while q is related to the sand flow in the wind direction. Here ignoring the correlated 
increase of these three parameters responding to the increase of wind force, we independently vary them as individual 
control parameters. 

As the initial condition of the numerical simulations, the height of 2D-CS crests are set uniform, i.e., hi{{)) — Ho{0 < 
i < N). By varying this initial height we control the amount of available sand in the field. Moreover, the initial wind 
directional position of crests are set wavy with small amplitude of sinuosity, that is, Xi{0) = HQ/20sm{2TTi/N). We 
check if the initial amplitude of sinuosity in Xi {0 < i < N) grows or not by measuring the quantity, 

N 

Var{t) =J2iMt) - x^+l{t))^ ■ 

i=l 

Simulations of the model exhibit three different phases depending on values of control parameter (Fig. 2]) : 

I) If Var{t)/Var{Q) < 10~^ is satisfied at t — 10* we consider the laterally extending straight transverse dune is 

stable. We call this 'ST-phase'. 

II) If Var{t)/Var{0) > 10"^ is satisfied at t = 10* we take the wavy shape of a transverse dune as temporally stable. 

We call this 'WT-phase'. 

III) If, at least, one 2D-CS shrinks to become hi{t) = within t < 10* thereafter the transverse dune will soon 
deform into the typical shape of barchan. We call this 'B-phase'. 




time 

FIG. 4: (a)~(d) Time evolutions from an initial straight transverse dune to a barchan. (e) Stable state of wavy transverse 
dune, (f) Typical time evolutions of the variance of dune crest line. Depending on the values of control parameters, three 
phases are observed; (I) ST-phase (real line): straight transverse dune is kept stable, (II) WT-phase (short dashed line): wavy 
transverse dune is kept stable, (III) B-phase (longer dashed line): initial transverse dune deforms into a barchan. 
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FIG. 5: Results of the numerical simulations varying two pairs of control parameters (a) Ho and q, (b) Vu and q. In white 
areas, ST-phase is realized, whereas, in grey areas WT-phase, in black areas B-phase is realized. Other parameters are set like 
(a) Vu = 0.1, Ud = 0.1, (b) ua = 0.1, Ho = 30.0. 
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FIG. 6: Results of the linear stability analysis of the fixed point of (|7a[) . (|7b| ) varying two pairs of control parameters: (a) 
{S,q), (b) {i/u,q)- White area indicate the conditions under which straightly extending transverse dune is stable, whereas 
in the black area such transverse dune is unstable. Other parameter are set like (a) ~ 5.0 x 10~^, Ud = 1.0 x 10~*, (b) 
= 1.0 X 10-"^, S = 15. 



In the first simulation, we fix j^„(= 0.1) and i'd{= 0.1) and vary two parameters: Hq and q (Fig. EJa)). The increase 
of Hq enhances the stabihty of the shape of straight transverse dune and the decrease of Hq destabilizes its shape to 
enforce the deformation to a barchan. This result qualitatively corresponds to the well known fact that barchans are 
formed in the field with small amount of available sand. 

In the next simulation, we fix = 0.1 and set Hq = 30.0 and vary two parameters, Vu and q (Fig. [5jb)). As f„ is 
set larger, the straight shape of transverse dune is more stabilized. On the other hand the increase of q destabilizes 
the transverse dune to enforce the deformation to barchan. In this way the balance between the wind directional flow 
and the lateral flow determine the stability of the shapes of transverse dunes. 

Finally, we conduct linear stability analysis of two 2D-CSs system which is considered as the simplest skeleton 
model of a transverse dune. Here the uniformly redistributing rule of escaped sand flow is applied like the above 
simulations. Thus, as long as two 2D-CSs are kept non-vanished, the total area of cross-sections S is kept constant. 
In this case, (pa|) and ([6b)) consist of four variables {hi, h2,xi,X2), and /i2 is expressed like /i2 — \/2AS — h\. Defining 
the wind directional relative crest position y = X2 — xi of two cross-sections, ([6a|) and (j6bp are simplified into two 
variables coupled equations. 



^ = (^(^^('^i) - 1) + ■^''(1-2) - ^.(1^2) + /r) (7a) 



Tt^^^y— —)^^^-^' 

-(Bi3d(i^2) - Ci3„(i^2)) (^^ + -1^ (7b) 

A fixed point {h\,y*) — (-\/A5',0) of (|7a |) . (ffbl) corresponds to the straightly extending state of a transverse dune. 
We analyze the stability of this fixed point. Figure. [SJa) and (b) show the linear stability of the fixed point in two 
sets of parameter spaces: (S', q) and (vu, q) spaces. In both figures, the fixed points are linearly stable in white areas, 
whereas they are unstable in black areas. These results qualitatively corresponds to those shown in Fig. [Sfa) and (b). 

In this letter, we proposed a skeleton model of 3D dune dynamics and studied the stability of the shape of transverse 
dunes through numerical and analytical methods. We also got a qualitative correspondence to the previous results 
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obtained through observations of real dunes and more comphcated simulation models. Beeause of the simplicity of 
the model, we expect this model supplies us with an effective tool for the theoretical study for the complex dynamics 
of 3D dunes. 
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